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Abstract 

Forests are important moderators of global atmospheric CO, emissions, making them a key 
focus of terrestrial C-cycling research. The 5th assessment report of the Intergovernmental 
Panel on Climate Change explicitly calls upon nations to enhance C-stock accounting and 
mitigate losses of global forest C sinks, which inherently will require more accurate and higher 
spatial resolution carbon accounting. Monitoring and predicting how disturbances, such as 
invasive species, influence forest C stocks and the transfer of C from live to dead pools 
remains a high priority both in the scientific and policy communities. We documented the 
effects of invasive North American beavers (Castor canadensis) on C-sequestration of riparian 
Nothofagus forests in Tierra del Fuego, Chile. Our paired plot sampling design quantified 
significant losses from beaver invasion in total aboveground, live standing, dead standing, and 
dead and downed C stocks (P< 0.001, paired ¢ tests). We extrapolated stand-level C losses to 
the entire study area using a Maxent habitat suitability model and estimated that 1.177 (+ 
0.103) Tg C would be lost if all of the predicted 18,384 ha of invasible habitat (14% of the total 
forested area) were colonized by beavers. These results document the impacts of invasive 
mammals on large terrestrial C sinks and highlight the need for understanding the magnitude 
of such effects across both landscape- and regional-scales. 


1 Introduction 


Forests play a pivotal role in global carbon cycling (Pan et al. 2011la; Smith et al. 2014), 
sequestering an estimated 25% of annual global anthropogenic C-emissions (Ciais et al. 2014); 
yet, the mechanisms underpinning forest C sink-source transitions remain poorly understood 
(Pan et al. 2011b; Running 2008). Assessing the relationships between forest disturbances and 
terrestrial C-cycling is paramount for improved carbon accounting (Grassi et al. 2017). 
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Disturbances including fires, pests, natural disasters, and pathogens can have profound effects 
on forest ecosystems and their long-term carbon capture and storage capacities, which often 
occurs through the conversion of forest carbon sinks to sources (Ellison et al. 2011; Trumbore 
et al. 2015; Turner 2010). Invasive species also act as major disturbance agents over large areas 
(Charles and Dukes 2008; Dukes and Mooney 2004) and are significant drivers of environ- 
mental change. Through the alteration of energy fluxes and biogeochemical cycles, forest- 
dwelling alien species pose a major threat to global forest stability (Tylianakis, 2008; Charles 
and Dukes 2008; Liao et al. 2008; Ehrenfeld 2010; Hicke et al. 2012). 

Forests comprise a recent major focus of international climate agreements due to their capacity to 
moderate atmospheric CO, concentrations (Berndes et al. 2016; Grassi et al. 2017; Houghton et al. 
2015). The 2015 Paris Agreement called for developed nations to begin making payments to 
developing nations as compensation for conserving forests as a part of REDD+ (Reducing 
Emissions from Deforestation, forest Degradation and other forest activities) and to report Nationally 
Determined Contributions (NDC) on a country-by-country basis (Grassi et al. 2017; IPCC 2008; 
UNFCCC 2017) that outline national mitigation objectives with forests as the basis for a quarter of 
emissions reductions through 2030. The effective implementation of REDD+ and NDCs hinges on 
more accurate and transparent measurement, reporting, and verification (MRV) of greenhouse gas 
inventories and C-flux from forests at the national level (Trumbore et al. 2015). Understanding the 
mechanisms responsible for forest C-flux provides the basis for predicting future atmospheric CO, 
levels, bolstering forest C sinks, and mitigating the effects of forest loss on C budgets (Trumbore 
et al. 2015; Pan et al. 201 1a, b). 

Invasive species can influence global biogeochemical cycling and can thereby represent a 
major global carbon emissions factor (Hicke et al. 2012). Nonnative forest pests such as the 
emerald ash borer, gypsy moth, mountain pine beetle, and hemlock wooly adelgid are causing 
massive C stock losses, as are invasive pathogens such as chestnut blight (Anagnostakis 1987), 
beech bark disease (Xu et al. 2016), and Dutch elm disease (Hicke et al. 2012; Loo 2009). 
While there is a growing body of literature on the effects of invasive insects and pathogens on 
forest C-sequestration, nonnative mammal impacts on forest C-cycling are largely lacking (but 
see Wardle et al. 2001, 2007). This signals the need for studies investigating such relationships, 
especially since the effects of invasive species are often unpredictable in their magnitude, 
direction, and duration. Reducing uncertainty about how nonnative biotic disturbance agents 
influence global forest C-cycling is a key step towards developing international agreements for 
monitoring forest C stocks and flux. 

In this study, we sought to develop an easily replicable protocol that could be used to rapidly 
estimate (at regional- and national-scales) the effects of an invasive mammal on aboveground forest 
C storage. We quantified the effects of the invasive North American beaver, Castor canadensis, 
on aboveground C stocks and changes to this C pool in response to land cover conversion from 
forest to wetland in riparian Nothofagus forests of Tierra del Fuego, South America. We also 
modeled the regional potential distribution of aboveground C pool transitions due to potential future 
beaver invasion. The study, from a remote and heavily forested region of the globe, adds to the 
growing body of literature on invasive species’ effects on forest C-cycling. 


2 Background 


North American beavers are ecosystem engineers that dramatically modify their habitat by 
clearing, damming, and flooding forested areas (Allen 1983). Within their native range, they 
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comprise a major riparian disturbance agent. As both an ecosystem engineer and herbivore, 
C. canadensis has the ability to affect both ecosystem function and species diversity (Nummi 
and Kuuluvainen 2013). Their impacts influence nearly all levels of the carbon and methane 
cycles within their native range, including streams, wetlands, and floodplains (Naiman et al. 
1988; Roulet et al. 1997; Whitfield et al. 2015; Wohl et al. 2012). 

Dramatic land cover alteration by beavers also triggers major ecosystem-level changes in 
this species’ invaded range in South America (Allen 1983; Anderson et al. 2009). The invasion 
of C. canadensis in South America has been recognized as the largest disturbance agent to the 
forests of Tierra Del Fuego (TDF) since the last ice age (Anderson et al. 2009). TDF is home to 
the southernmost forest on the planet, and it was recently named one of the most pristine 
wilderness areas in the world because of its extensive size, remoteness, and lack of human 
disturbance (Mittermeier et al. 2002). 

Beavers were introduced to Isla Grande, TDF in 1946 by the Argentinian government for 
the fur trade industry (Anderson et al. 2009). They have since expanded to most of the islands 
in the archipelago and, more recently, to mainland South America (Pietrek et al. 2017; 
Westbrook et al. 2017). Their population has expanded exponentially in the absence of 
predators and with abundant available habitat (Lizarralde 1993; Westbrook et al. 2017). 
Current C. canadensis population densities in TDF range from 0.2 to 5.8 colonies kny? 
compared to a range of 0.08 to 1.4 colonies km~ within their native range (Anderson et al. 
2006; Gurnell 1998; Lizarralde 1993; Westbrook et al. 2017). While their dam building can 
have myriad negative impacts on Tierra del Fuego riparian ecosystem function, including 
shifts in understory succession (Pastur et al. 2006), pond formation resulting in methane 
emissions (Whitfield et al. 2015), and changes in inorganic and organic soil nitrogen concen- 
trations (Lizarralde et al. 1996), neither the local- or regional-scale effects of beaver invasion 
on forest carbon sequestration have been estimated in South America. 


3 Methods 
3.1 Study area description 


The study was conducted in Karukinka Park (53° 42’ 0" S, 69° 19’ 0” W), which is owned and 
managed by the Wildlife Conservation Society. Karukinka is located on the Chilean side of 
Tierra Del Fuego’s main island, Isla Grande, and is approximately 300,000 ha in size (Tapia 
2010) (Fig. 1). The study area consists mainly of sphagnum bogs and beech forests, both of 
which play an important role in carbon capture and storage (Saavedra et al. 2011). The forests 
of Karukinka are dominated by two deciduous tree species, Nothofagus pumilio and N. 
antarctica, and one evergreen species, N. betuloides which dominates the coastline (Barrera 
et al. 2000; Pisano 1977; Saavedra et al. 2011°). Mean annual temperatures range from 2.6 to 
6.6 °C, with mean monthly temperatures of 9.5 °C in January and 0 °C in July (Hobbs et al. 
1998; Tuhkanen 1992). The average annual precipitation ranges from 408 to 620 mm (Pisano 
1932): 


3.2 Species biology 


The native range of C. canadensis includes all of North America, except for the arctic tundra 
and Florida (Jenkins and Busher 1979). Beavers are generalized herbivores, but they show 
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Fig. 1 Map of the study area and forest cover (Tapia 2010) within Karukinka Park, WCS, Tierra del Fuego, Chile 


preference for particular tree species in both their native range in North America and their 
invaded range in South America. In Tierra del Fuego, they have a preference for N. pumilio 
and avoid Drimys winteri and Maytenus magellanica (Wallem et al. 2007). Bank slope, stream 
gradient, hardwood canopy cover, and stream width and depth are important factors influenc- 
ing beaver dam construction and subsequent colony formation for C. canadensis in their native 
range (Howard and Larson 1985; McComb et al. 1990; Slough and Sadleir 1977). However, 
the extent to which these factors influence beaver dam construction and colonization within 
this species’ invaded range is not well documented. 


3.3 Beaver presence data 


Beaver dam presence data were generated using Google Earth Pro satellite imagery (n = 2068 dams) 
for Karukinka Park. Google Earth Pro imagery was used exclusively to identify points where 
beavers were present for subsequent maximum entropy (MaxEnt) (Phillips et al. 2006) species 
distribution modeling and to identify dams to be located for field data collection. Dams were easily 
identified in the imagery by the characteristic triangular shape of water behind each dam which was 
surrounded by an irregular perimeter of standing dead trees that had died from intermittent dam 
flooding events during the period of beaver dam occupation. The match between beaver presence in 
the field and the aerial triangular shape of the dam was also validated via field visits and drone flights 
by Karukinka Park manager, Daniel Teran, in 2015 prior to satellite imagery interpretation (Fig. 2). 
Dams were defined as single barriers to water flow created by C. canadensis. Each individual dam 
was marked at its centroid in the imagery. We used cloud-free satellite imagery for the study area 
from Google Earth Pro for 2004, 2014, and 2015. Stands of dead riparian trees resulting from other 
causes were also evident in the imagery. However, such sites lacked the characteristic triangular 
beaver dam shape and were thereby excluded in the generation of the beaver dam presence dataset. 


3.4 Raster environmental data 


A suite of raster environmental data layers was derived from 30-m-resolution ASTER GDEM V2 
(USGS NASA Earth Explorer) digital elevation data for the spatial extent of Karukinka Park. 
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Fig. 2 Aerial imagery of dams created by C. canadensis in Tierra del Fuego, Chile depicted in a aerial image of 
beaver dam taken by a drone (Daniel Teran, unpublished data), b the same beaver dam in Google Earth Pro 
imagery, and c other beaver dams identified in Google Earth imagery (Google Earth Pro, Digital Globe 2017) for 
the beaver dam presence points in the MaxEnt modeling procedure. Standing dead trees surrounding the beaver 
dams themselves occurred in response to intermittent flooding events during the period of beaver dam occupation 
and were therefore attributable to beaver activity 


Elevation, upstream flow, downstream flow, distance to stream, and slope were all produced in 
ArcMap 10.4.1 (ESRI 2017) and were projected to WGS 1984 UTM 19S (Table 1). The raster data 
layers included in the modeling procedure were chosen based upon their relevance to riparian 
habitats and upon the previously mentioned factors that influence beaver distributions within their 
native range (Howard and Larson 1985; McComb et al. 1990; Slough and Sadleir 1977). A 
Landsat 8, 30-m-resolution vegetation map by Tapia (2010) was also used as an environmental 
predictor variable with six vegetation cover-types: bog, closed forest, open forest, meadow and 
pasture, rocky outcrop, and bodies of water. We were conservative in the environmental predictors 
we used in our model; we included only biologically meaningful variables and that would be 
indicative of C. canadensis activity in the model construction. 


3.5 Habitat suitability modeling 


We utilized the maximum entropy (MaxEnt) algorithm (Elith et al. 2011; Phillips et al. 2006), 
version 3.3.3k to produce a habitat suitability map for C. canadensis within the study area to assess 
the invasion potential of this species in Karukinka Park. MaxEnt is a well-known and high- 
performing niche modeling method (Elith et al. 2011; Jarnevich and Reynolds 2011; Ward 2007) 
that uses a simple deterministic algorithm to identify the most suitable habitat across a given area 
with specified environmental covariates. MaxEnt has previously performed well in modeling 
species distributions of aquatic invasive species, especially since it is a presence-only niche modeling 
algorithm (Bouska et al. 2015; Hayer et al. 2014; Leidenberger et al. 2015; Poulos et al. 2012). We 
acknowledge that one limitation of invasive species niche models is that species on the invasion 
front comprise the species’ realized niche, but that the species may not yet occupy their fundamental 
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Table 1 Environmental predictors with their native spatial resolution and source. All predictors were included in 
the habitat suitability model 


Layer Type Description Spatial Source 
resolution (m) 


Elevation Raster-continuous Acquired from NASA USGS in 30 ASTER GDEM 


Earth Explorer USGS! 
Upstream flow Raster-continuous Upslope distance flow path in the 30 ASTER GDEM 
Flow Length toolbox in ArcGIS USGS 
Downstream — Raster-continuous Downslope flowpath in the Flow 30 ASTER GDEM 
flow Length toolbox in ArcGIS USGS 
Distance to Raster-continuous Euclidean distance to flow accumulation 30 ASTER DEM 
stream raster from ArcHydro extension USGS 
of ArcGIS 
Vegetation Raster-categorical Six distinct vegetation types as 30 Tapia 2010 
derived by Tapia (2010) 
Slope Raster-continuous Surface toolbox in ArcGIS 30 ASTER GDEM 
USGS 


' ASTER GDEM United States Geologic Survey: ASTER GDEM is a product of METI and NASA. Accessed 
through: https://earthexplorer.usgs.gov 


niche space within the invaded habitat (e.g., Elith and Leathwick 2009; Rédder and Lotters 2009). 
We accounted for this limitation in this study in our choice ofa presence only algorithm and through 
the incorporation of a bias file into the modeling procedure, as described below. MaxEnt uses 
presence only data and takes random background points as pseudo-absence. The use of background 
points as pseudo-absence could negatively affect the model performance due to selection bias. 
Selection bias can happen in different ways; it can come from presence data that is biased in a study 
area because of ease of accessibility to certain areas (because of roads or because an area is near to 
shelter) or through the placement of pseudo-absence points in areas that have not yet but could be 
invaded by beavers. Our use of Google Earth Imagery for identifying beaver dam presence points 
allowed accessibility to remote areas and thus reduced biases associated with access to remote areas 
of the Karukinka Park. However, we also attempted to reduce biases associated with the placement 
of pseudo-absence points in potentially invasible areas by generating a bias file in ArcMap 10.2 
(ESRI 2017) that limited the extents of pseudo-absence points to areas that were > 200 m from a 
beaver dam presence location, thereby making pseudo-absences a better approximation of actual 
absences. 

Maxent evaluates model performance by calculating area under the curve (AUC) of the 
receiver operator characteristic plot. The AUC is a threshold-independent measure of model 
performance and ranges in value from 0 to | with 1 being the highest performance possible. 
Values greater than 0.9 indicate high accuracy, values within the range 0.7—0.9 indicate good 
accuracy, and values below 0.7 indicate low accuracy (Swets 1988). 

The final model was run 25 times with random substitution and replacement of training and 
testing data (55% training and 45% testing subsampling). Maps of the mean performance and 
standard deviation from the mean over the 25 runs were produced following Jarnevich and 
Reynolds (2011). We applied a threshold to the final MaxEnt habitat suitability map using 
equal training sensitivity and specificity to create a binary presence/absence map, which is 
considered a robust thresholding approach (Liu et al. 2005). We used this binary map to 
determine the total area vulnerable to beaver invasion in hectares within the study area. We 
further categorized this area by vegetation type and used the total hectares of vulnerable forest 
to estimate potential aboveground C loss by beaver invasion. 
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3.6 Quantifying carbon loss 


Riparian dams were first identified via Google Earth imagery. Dams for field sampling were chosen 
based upon accessibility and whether they were still occupied and therefore flooded. Unoccupied. 
dams were selected for field measurements so that surface coarse woody biomass could be 
sampled. Dams were located in the field using a handheld GPS unit, and plot centers were chosen 
randomly upon arrival by picking a random azimuth from the center of the dam and walking ten 
paces in that direction. We employed a paired plot sampling design where plots were placed at an 
invaded dam site and then in adjacent, intact, and un-invaded forest stand. Adjacent, intact, 
forest stands were selected by walking more than 50 m into the intact forest away from the 
beaver dam. Edge effects of the beaver dam on forest structure were minimized by ensuring that 
plots were placed well within the living closed canopy forest so that no portion of the intact 
forest plot was influenced by the nearby beaver dam. 

Circular, variable-radius plots were used to inventory all live and standing dead trees using a 
metric, two-basal area factor prism. The only species encountered in our study was N. pumilio, and 
thus, standing C stocks were estimated for each plot using an allometric equation for this species 
according to Schmidt et al. (2009). Plot radii ranged from 10 to 18 m depending on forest density so 
that intact forest plots contained between 40 and 50 trees each. For each tree in each plot, we 
recorded diameter at breast height (DBH), species identity, and condition (standing dead or living). 
Coarse woody debris, or surface C stocks, were estimated using point relascope sampling with a 
20° angle of inclusion according to Gove et al. (1999). The diameter of each log end and its length 
was recorded for each tree that fell within the inclusion angle (Gove et al. 1999). Biomass was 
calculated following Gove et al. (1999), and surface C from coarse woody debris was estimated 
using the bulk density of N. pumilio reported in Schmidt et al. (2009). Carbon stocks were 
estimated for standing live trees, standing dead trees, and dead and downed trees, and total carbon 
was calculated in megagrams per hectare in all paired (invaded and un-invaded) plots. Differences 
between C stocks in invaded and un-invaded forest stands were statistically assessed via paired t 
tests in R (R Development Core Team 2018). 

Aboveground C-storage change due to C. canadensis invasion was estimated by first calculat- 
ing the mean difference in C stocks between the paired (7 = 34 pairs) un-invaded and invaded 
stands. We estimated C loss for the entire study area by multiplying the total vulnerable forested 
area (from the Maxent model in hectares) by the difference in C stocks (Mg ha“') as an estimate of 
C losses if the entire potential habitat was invaded by beavers. As a more conservative estimate of 
loses, we also estimated C loss if beavers were to invade just half of the potential invasible habitat. 


4 Results and discussion 
4.1 Beaver habitat suitability modeling 


The MaxEnt model performance revealed a mean AUC (area under the curve of the receiver 
characteristic plot) value across all runs was 0.846 (SD = 0.004), indicating good overall model 
performance predictors including distance to a stream, vegetation cover, and slope as the most 
important predictors of C. canadensis habitat, consistent with its known habitat requirements. 
Other minor contributors to the model (> 10%) included upstream flow, elevation, and downstream 
flow (Fig. 3; Table 3). The Maxent environmental response curves also identified sites in 
Karukinka Park that were suitable for beaver dam establishment (Figs. 4 and 5). These 
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characteristics were consistent with known habitat selection preferences in the native and invaded 
range (Anderson et al. 2009; Howard and Larson 1985; Lizarralde 2004; McComb et al. 1990; 
Slough and Sadleir 1977) suggesting that flat, lowland habitats near streams are a primary forest 
conservation and beaver eradication target. In addition to providing an estimate of the total 
invasible area, the habitat suitability model for Karukinka can also be used as a decision support 
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Fig. 3. Environmental response curves for all environmental predictor variables included in the final Maxent 
potential C. canadensis distribution model for Karukinka Park, Chile. Red lines indicate mean values for the 25 
iterations of the model. Blue shading indicates + SD. Probability by vegetation type is shown in the bottom left 
panel in red + one (dark blue) — one (light blue) SD 
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Fig. 4 Maxent potential distribution maps for C. canadensis in Karukinka Park, Chile. a The habitat prediction 
map derived from the mean of the 25 iterations. b The standard deviation among the 25 iterations using different 
subsets of the presence points for training (55%) and testing (45%). ¢ The presence—absence map after 
thresholding. Maximum sensitivity plus specificity was used as the threshold to derive the binary prediction 
map product 


tool for informing efforts focused on limiting the spread of beavers and for guiding eradication 
efforts at the invasion front (Elith et al. 2011). 
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Fig. 5 Forests (Tapia 2010) that are vulnerable (18,384.48 ha) to invasion by C. canadensis in Karukinka Park, 
Chile as predicted by Maxent 


4.2 Effects of beaver activity on aboveground C stocks 


North American beaver activity in Tierra del Fuego has triggered significant losses of Nothofagus 
aboveground C stocks via land cover conversion from forest to wetland. The aboveground field 
carbon inventories revealed that beaver-invaded areas experienced significantly lower total above- 
ground C stocks relative to adjacent intact forests (paired ¢ test, df= 33; fvalue = 208.8; p value < 

0.0001). Standing dead (df= 33; f-value = 11.65; p value = 0.0017), dead and downed (df= 33; f 
value = 34.43; p value <0.0001), and total carbon (df= 33; f value = 137.97; p value < 0.0001) 
stocks also differed significantly between invaded and uninvaded areas (Table 2). Mean C conver- 
sion from live to dead pools due to C. canadensis invasion was 64+5.6 Mg ha, which we 
estimated by taking the difference in C stocks between the paired invaded and uninvaded forest 
stands. 


Table 2 Effect on aboveground carbon from invasion of C. canadensis in riparian forests. Mean (+ SE) 
aboveground carbon storage in both invaded and un-invaded riparian forests were estimated in un-invaded 
(n = 34) and invaded (n = 34) plots in Karukinka Park, Tierra del Fuego, Chile 


Un-invaded (Mg ha™!) Invaded (Mg ha™!) 
Living aboveground 80.7 (£5.43) 2.01 (+ 1.45) 
Standing dead 9.65 (£1.99) 20.96 (+ 2.68) 
Coarse woody debris 2.68 (£0.31) 5.83 (£0.48) 
Total carbon 93.05 (£5.69) 28.81 (£2.82) 
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Our results suggest that beaver invasion can result in major differences between aboveground C 
stocks in invaded versus un-invaded forest stands. However, they also highlight differences in the 
compartmentalization of the aboveground C stock. In un-invaded stands, most of the C was stored 
in live standing trees, while in beaver-invaded and recently abandoned beaver colonies, the 
majority of the C stock was stored in standing dead and dead and downed trees. Other authors 
have observed that post-abandonment beaver ponds do not revert back to Nothofagus forest stands 
over the long-term, but rather convert to sedge- and grass-dominated meadows that persist 
indefinitely (Anderson et al. 2009; Pastur et al. 2006). While the soils of the plant communities 
that form after dam abandonment and pond drainage are C-rich following the incorporation of 
coarse and fine organic matter from the forest die-off event, they also comprise significant carbon 
dioxide and methane greenhouse gas emissions “hot spots” (Scanlon and Moore 2000; Updegraff 
et al. 1995; Vecherskty et al. 2011), or emissions factors. 

The prolonged inundation of forest stands from beaver damming and subsequent sedimen- 
tation and nutrient deposition amplifies these effects and prevents tree regeneration from the 
seed bank following beaver colonization (Anderson et al. 2009; Pastur et al. 2006). Although 
other major forest die-off events release carbon slowly from dead trees, (Moore et al. 2013), 
the rapid conversion of Nothofagus forests in TDF from intact forest to flooded beaver dam 
precludes rapid forest recolonization and recovery at these sites, at least over the short-term. 
Our field data and observations indicate that tree recruitment has yet to occur at older, 
abandoned dam sites, even 10 years after beaver dam abandonment. This suggests that beaver 
invasion stimulates a persistent loss of riparian forests and that future invasion of beavers will 
likely have major consequences for forest C sequestration across the region. 

While this study documents the conversion of an aboveground C from a sink to a source, we also 
recognize that this process likely occurs over time-scales spanning from months to decades. 
Flooding caused by dam creation kills trees over short time-scales (i.e., months), but the loss of 
carbon and CO), emissions to the atmosphere can occur over longer time-scales, depending on the 
permanency of the dam (Johnston 2014, 2017). Standing dead snags can persist for years after tree 
death, and the pool created by the beaver dam subsequently traps and submerges woody debris and 
sediments in an anoxic environment that can store organic C for years to decades post-dam 
abandonment (Johnston 2014; Naiman et al. 1988, 1994). The process by which stored C-stocks 
in standing live trees are either released into the atmosphere as CO) after tree death or remain trapped 
in underwater organic sediments is complex. Our knowledge of how and when such changes occur 
is largely restricted to northern latitudes within the native range of the species. Such processes 
warrant further study within this species’ invaded range as a means of elucidating the long-term 
consequences beaver invasion on C-cycling. Longitudinal studies that examine how long beaver 
dam sites remain ponds that can trap organic sediments after abandonment could facilitate our 
understanding of the rate of C sink-source transitions, as could research that evaluates pond— 
meadow-—pond transitions in response to the recolonization of former dam sites by beavers. 


4.3 Landscape-scale aboveground C losses 


Our Maxent habitat suitability model predicted that 18,384 ha (14% of the total forested area) of 
forested land is vulnerable to invasion by C. canadensis (Fig. 4; Table 3). If the entire potential 
invasible forested area (18,384 ha) of Karukinka was invaded by beavers, 1.177 (+0.103) Tg C, or 
10% of the total carbon stored in the forests of the park, would be lost based on the potential 
invasible area and our field-based estimates of C storage). A more modest calculation of 0.588 (+ 
0.051) Tg C would be lost if only half of the potential invasible area were colonized by 
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Table 3 Relative contribution of the environmental covariates to the Maxent model. Percent contribution reports 
the gain by including that specific variable at each step of the Maxent algorithm. The permutation importance 
reports the percent contribution of each environmental covariate 


Environmental predictor % Contribution Permutation importance 
Distance to stream 42.1 35 

Vegetation 24 20.1 

Slope 13.9 19.1 

Upstream flow 9.4 12 

Elevation 5.4 6.9 

Downstream flow 23 6.9 


C. canadensis. We acknowledge that it is unlikely that beavers would invade all of the predicted 
suitable habitat, but the invasion of even half of the potential habitat would result in significant forest 
conversion and aboveground C stock losses. To put these losses into context, between 2006 and 
2010, C loss due to deforestation for the entire conterminous United States (which contains 
approximately 1800 times the forested area as is within Park Karukinka) was 6+1 Tg C year"! 
(Harris et al. 2016), highlighting the magnitude of this potential disturbance on global C- 
sequestration due to beaver activity. 


5 Conclusion 


Our study presents a novel method for quantifying the effects of invasive species on terrestrial 
C-cycling by combining a habitat suitability model with field-based C stock data—an ap- 
proach that could be easily applied elsewhere to quantifying the effects of other invasive 
species on aboveground C budgets. Modeling the species’ range provided a means for 
estimating the potential large-scale effects of beaver invasion on forest carbon sequestration 
and allowed us to extrapolate field C-inventory data across a wide geographical area. Such an 
approach provides a mechanism for generating landscape- and regional-scale estimates of the 
conversion of C from live to dead pools in response to mammal invasion, which is useful for 
newly-mandated national carbon-accounting in this region and elsewhere. 

The study also highlights the magnitude of change that an invasive mammal can have on forest 
C-cycling and storage. The estimated loss of above-ground forest C caused by beaver invasion, 
combined with poor post-invasion forest regeneration potential (Pastur et al. 2006), suggests that 
active restoration of previously invaded areas will be necessary for mitigating the effects of beaver 
invasion on South American C-storage. While this study estimates aboveground C loss, further 
research could contribute to a greater understanding the effects of beavers, as well as other 
invasive mammals on aboveground C-storage and biogeochemical cycling within the ecosystem, 
including in-stream effects and effects on the archipelagos peat bogs that result from forest 
conversion by beavers. Globally, including South America, C. canadensis has already been 
recognized as a significant source of methane emissions via pond formation (Whitfield et al. 
2015), although data are lacking from within this species’ invaded range. Thus, a combined 
mitigation approach that focuses on beaver spread prevention, population reduction through 
beaver eradication efforts, and forest restoration will be necessary for ameliorating the negative 
effects of invasive beaver activities on C budgets (Anderson et al. 2009; Pastur et al. 2006). 

Preventing and reducing deforestation has been specifically identified by the IPCC as one of the 
largest and most immediate ways to decrease global CO, emissions in the short term (IPCC 2008; 
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Smith et al. 2014). The signing of the Paris Agreement has placed increased international attention 
on the importance of forests in reducing atmospheric CO . The successful implementation of 
REDD+ and reduction of emissions by a quarter through 2030 depends on more accurate and 
transparent C-monitoring, reporting, and verification. More detailed and process-based studies of 
the mechanisms responsible for alterations to the C cycle within forests are needed for generating 
accurate forest C stocks and emissions factors. Not only does C. canadensis pose a serious threat to 
the forests of TDF but they also threaten mainland South American forests, and the effects of the 
widespread invasion of this invasive mammal could have effects on C-cycling that spread well 
beyond the study area. Thus, incorporating invasive species as a C emissions factor in land use/land 
cover change carbon accounting represents an important step in achieving more accurate IPCC Tier 
II and Tier III C-accounting standards. 
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